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BLAST: A FAR-INFRARED MEASUREMENT OF THE HISTORY OF STAR FORMATION 
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ABSTRACT 

We directly measure redshift evolution in the mean physical properties (far-infrared luminosity, 
temperature, and mass) of the galaxies that produce the cosmic infrared background (CIB), using 
measurements from the Balloon-borne Large Aperture Sub- millimeter Telescope (BLAST), and Spitzer 
which constrain the CIB emission peak. This sample is known to produce a surface brightness in 
the BLAST bands consistent with the full CIB, and photometric redshifts are identified for all of the 
objects. We find that most of the 70 /im background is generated at z < 1 and the 500 /im background 
generated at z > 1. A significant growth is observed in the mean luminosity from ~ 10 9 -10 12 L Q , 
and in the mean temperature by 10 K, from redshifts < z < 3. However, there is only weak 
positive evolution in the comoving dust mass in these galaxies across the same redshift range. We 
also measure the evolution of the far-infrared luminosity density, and the star-formation rate history 
for these objects, finding good agreement with other infrared studies up to z ~ 1, exceeding the 
contribution attributed to optically-selected galaxies. 

Subject headings: cosmology: observations — cosmology: diffuse radiation submillimctcr — galax- 
ies: evolution — galaxies: starburst 



1. INTRODUCTION 

The spectrum of the diffuse extragalactic background 
radiation reflects the physical processes that have dom- 
inated the evolution of structure in the Universe. The 
Cosmic Microwave Background, the relic radiation of the 
early Universe, is the dominant contribution. 

The second most important component consists of 
two broad peaks at around 200 /im and 1 /im which 
carry roughly equal energy, and are presumably asso- 
ciated with light emitted by stars throughout cosmic 
time. In the far- infrared (FIR), the isotropic Cosmic In- 
frared Background (CIB) was first detected with the FI- 
RAS (IPuget et alj[l99l iFixsen et al][l99l and DIRBE 
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(jHauser et al.l I1998D instruments on board the COBE 
satellite. It consists of thermal emission from warm dust 
which enshrouds star-forming regions in galaxies; thus, 
a detailed account of the CIB will by necessity include 
constraints on the history of star formation in the Un i- 
verse (see also lDwek eFal]H99llHauser fc Dwekll200lD . 
In strong contrast to the emission at w 1 /jm, the warm 
dust is optically thin, so its intensity reveals the entire 
mass without assumptions or corrections. In this paper 
we resolve the contribution to the CIB as a function of 
redshift. This constitutes a measurement of the FIR his- 
tory of the Universe which is closely related to the history 
of star and galaxy formation. 

Large fractions of the CIB have been resolved into con- 
tributio ns from individual so urces at 24 /im, 850 /im, and 
1.1mm ([La gachc et al. 2005, and references therein), but 
the background is several tens of times smaller at these 
wavelengths than it is at its peak. Near the peak, cur- 
rent and anticipated FIR and sub millimeter experime nts, 
including the SPIRE instrument (jGriffin et al.ll2003h on 
the Herschel satellite, cannot resolve the CIB into in- 
dividually detected sources because of confusion arising 
from the finite instrumental angular resolution. 

It is, however, possible to study the average contri- 
bution that a given class of objects makes to the back- 
ground with a stacking analysis consisting of calculating 
the covariance of known source positions with maps from 
these experiments. This ap proach has been successfully 
used by many authors (e. g. iDole et all 12006: Dy e et all 
120071: iSerieant et all 12008ft . Individual sources brighter 
than 60/iJy detected in deep surveys with the Multi- 
band Imaging Photometer for Spitzer ( MIPS) at 24 
resolv e most of the CIB at 24 /im (~ 70%. lPapovich et alJ 
2004). Stacking analyses have shown that MIPS sources 



2 



Pascale et al. 



consitute a larg e fraction of the background at 70 and 
160 /an (~ 80%, iDole et alf2 006). The bulk of this emis- 
sion is from comparatively low redshift galaxies, while at 
850 /im, and 1 mm the CIB is dominated by massive star- 
forming galaxies at higher redshifts (z « 2.5). Besides, 
iDve et all (|2007h determine that these extreme galaxies 
contribute relatively little emission to the important re- 
gion near the CIB peak which forms the transition be- 
tween plentiful nearby galaxies and rarer distant sources 
of the submillimetre regime. 

BLAST has observed an 8.7 deg 2 region encompass- 
ing the Extended Chandra Deep Field South (ECDFS) 
and Great Observatories Origins Deep Survey South 
(GOODS-S). Additional time was spent observing a 
smaller 0.9 deg 2 area at the center of this field, though the 
coverage is still broader than most of the other deep data 
available at other wavelengths. We refer to these fields 
as BLAST GOODS-S Wide (BGS-Wide) and BLAST 
GOODS-S Deep (BGS-Deep) respectively. We note that 
the variance of the BGS-Deep maps are dominated by 
the confusion from point sources arising from the 36"- 
60" BLAST beams, rather than instrumental noise. 

Observing in three broad bands centered at 250, 350, 
and 500 /im, BLAST provides a new data set of unique 
size, resolution, and spectral coverage to study the prop- 
erties of the CIB on the submillimeter side of its maxi- 
mum. Together with MIPS observations at 70 /im, it is 
now possible to constrain the CIB without assumptions 
about any of the physical properties of the sources con- 
tributing to it. 

The possibility of determining number counts, spectral 
energy distributions, and clustering properties of sub- 
millimetre sources directly from correlations within the 
maps, without the requirement to first ex tract a cat- 
alog o f point sources, was pointed out by iKnox et al.l 
(2001). The approach has the important advantage of 
avoiding flux and number densi ty biases associated with 
the process of catalog s election. iDevlin et al.l (|2009t ) and 
iPatanchon et al.l (|2009l ) use this approach to estimate the 
number counts of extragalactic sources in the BLAST 
bands acr oss nearly fi v e orde rs of magnitude in flux den- 
sity, and iViero et ail j2009) measure the clustering of 
B LAST source s via th e power spectrum 

IDevlin et all (|2009l ). iMarsden etafl (|2009h . and this 
work extend that idea by measuring the covariance with 
external data sets, but again with out forming a B LAST 
catalog, which is insted studied bv lDve et alj (|2009f ) who 
find MIR and radio counterparts for a sample of bright 
BLAST sources and measure their redshift distribution. 
Using this stacking technique, it is found that the to- 
tal intensity emitted by MIPS 24 /im-selected sources at 
the BLAST wavelengths is compatible with the FIRAS 
measurements of the CIB intensity within the FIRAS 
la uncertainty, which is w 25%. The 24 /im-selected 
source redshift distribution has a median of z ~ 0.9 (see 
$3]), and it is certainly possible that a population with 
fainter 24 /im flux densities is missing, possibly associated 
with the submillimeter galaxies detected by SCUBA and 
MAMBO. 

We extend this work by finding the covariance with 
catalogs which contain redshift information. We present 
a detailed analysis of all sources that contribute signifi- 
cantly to the CIB at 70 /im, 250 /im, 350 /im, and 500 /im 
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Fig. 1. — Extent of the catalogs in the BGS-Deep region ob- 
served with BLAST. The area used for the stacking analysis is the 
intersection of the FIDEL and MUSYC regions and extends over 
590arcmin . 

as a function of redshift, even those well below the confu- 
sion limit. From this we derive the evolutionary history 
of the CIB. The results from these five analysis papers 
will be used by Chapin et al. (in preparation) to infer 
the evolution of the FIR luminosity function. 

A portion of the BGS-Deep region has been surveyed at 
many wavelengths from which we have derived an almost 
complete set of redshifts (~ 95%) for MIPS sources. By 
exploiting BLAST and MIPS 70 /tm data, we estimate 
the mean physical parameters for these galaxies (FIR lu- 
minosities, temperatures, and dust masses), and the total 
FIR luminosity evolution and star formation history of 
the Universe. Our results are in good agreement with ex- 
isting measurements. We find that stars form in optically 
obscured dusty galaxies at rates ~ 3 times larger than 
those previously estimated from optical data at redshifts 
up to z ~ 1. We also estimate the variation of the total 
comoving dust mass density with redshift, finding only 
moderate evolution, within the experimental uncertain- 
ties, over the range z < 3. 

A flat cosmological model with £l\ — 0.7, h — 0.7 has 
been used throughout. 

2. DATA 

In this section we describe each of the data sets used. 
Their locations and angular extents are depicted in Fig- 
ure Q] 

2.1. BLAST Submillimeter Imaging 

The BGS-Deep region observed with BLAST covers 
the central 0.9 deg 2 of a shallower map imaging the full 
8.7 deg 2 area in BGS-Wide. 

The raw time ordered data (TOD), consisting of 
a set of voltage time-streams fr om each of the 
BLAST detectors, are pre-processed (|Truch et al.l 120081 
IPatanchon et al.ll2008l ). Any corrupt samples are flagged, 
and the data are deconvolved using the instrumental 
transfer functions. The TODs are b inned into maps us- 
ing the telescope pointing solution (|Pascale et al.H200"8f) 
whi ch has an estimated R MS positional error of less than 
5" (jMarsden et al.ll2008D . The absolute photometric cal- 
ibration has an estimated accuracy of ~ 10%, which is 
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highly correlated across the BLAST bands (|Truch et al.l 
2009). Calibration errors and color corrections are ac- 
counted for in our SED fits. 

The maps used in this work have been generated using 
OPTBIN (Pascale et al. in preparation). The algorithm 
performs common-mode suppression and high-pass fil- 
ters the output above the low-frequency "knee" of the 
detector noise. The filtered data streams are then pro- 
jected onto a map. This treatment is appropriate for 
point source extraction as it does not distort the data 
on timescales comparable to the beam-crossing times. 
OPTBIN reduces the 100 hour long BGS-Wide data-set at 
250 /im in one hour on a single-CPU machine, and it is 
also used as part of an end-to-end instrumental simulator 
which allows Monte Carlo analysis. 

Results have been compared and found to be com- 
patible with those obtained using maps reduced with 
the more comput a tiona lly intensive SANEPIC algorithm 
(|Patanchon et all I2008D . which provides the least- 
squares solution to the map-making equation. 

2.2. Spitzer Data 

Extensive Spitzer maps cover most of the BGS-Deep 
and Wide regions. The Spit zer Wide- Area Infrare d Ex- 
tragalactic survey (SWIRE, lLonsdale et all [20031 . data 
release 3), provides ^7.5 deg 2 maps well matched to the 
BGS-Wide coverage. We use the BLAST data in con- 
junction with the 70 fim SWIRE data to constrain the 
brightness of the peak in the CIB. 

Deep imaging of this field is available in the four 
IRAC (Infrared Array Camera) bands from the SIM- 
PLE (Sp itzer's IRAC and M USYC Public Legacy of the 
ECDFS, iDamen et all I2009D survey. MIPS data from 
the Spitzer Far-Infrared Deep Extragalactic Legacy (FI- 
DEL) survey are also available. The ca talog of 24 /im 
sources detected in the FIDE L maps of Ma.gnclli et alJ 
(l2009h is the s ame a s used by iDevlin et al.1 (|2009l ) and 
Marsdc n~et~al] (|2009D . and utilizes the IRAC sources de- 
tected in SIMPLE as a positional prior. 

The faintest 24 /im source in the catalog has a flux 
density of 13/iJy, but the detections at S24 < 30/iJy 
should be considered tentative. The catalog includes 
9,110 sources covering 700 arcmin 2 , with IRAC flux den- 
sities (at 3.6, 4.5, 5.8, and 8/im) available for each entry. 

2.3. Optical Photometric Redshifts 

We use five different catalogs of photometric redshifts 
available in the BGS-Deep region to assign a redshift to 
t he majority of the MIR sources in the FIDEL catalog. 

Bramme r et al.1 (|2008f) released two catalogs of photo- 
metric redshifts in the BGS-Deep region u sing MUSYC 
(MU ltiwavelength Survey by Yale-Chile. | Tavlor et al.l 
l2009h and FIREWORKS (jWuvts et al.1 I2008H imaging 
data, calculated with a new algorithm called eazy. eazy 
performs linear combinations of SED templates to find 
the best fits to measured photometry, which in these two 
cases span the optical (UBVRI) and the IR (JHK), as 
well as the 4 IRAC bands. They provide values for the 
Q z statistic which indicates the confidence level associ- 
ated with the inferred redshifts. As explained therein, 
redshifts with Q z < 3 result in a statistical error of 
<t 2 = Az/(1 + z) < 0.1. Whenever the Q z statis- 
tic is available, we select redshifts using this cut. The 
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Fig. 2. — The redshift distribution of MIR sources is shown as a 
solid line for entries having a reliable photometric redshift assigned 
to them. The distribution peaks between 0.5 < z < 1 and exhibits 
a substantial tail up to z ~ 2. The dashed histogram shows the 
distribution of sources with redshifts derived from the IRAC flux 
density ratios, which tend to cluster at 2 > 1. 

MUSYC catalog covers an area of 900 arcmin 2 and con- 
tains - 17,000 sources, while FIREWORKS extends over 
a 140 arcmin 2 area and has ~ 6,300 sources. The typical 
redshift errors are a z < 0.04. 

The GOODS-MUSIC (|Grazian et al.l l200l catalog 
contains ~ 15,000 sources in a region overlapping FIRE- 
WORKS. About 6.5% of the redshifts listed in the cat- 
alog are spectroscopic, and the others are photometric 
estimates. The photometric redshifts ar e based on data 
at sim ilar wavelengths to those used by iBrammer et al.1 
(120081). but the algorithm itself differs from eazy. 
I Taylor et alJ (|2009D show that the two sets of redshift es- 
timates are in equally good agreement with the available 
spectroscopic redshifts. However, this agreement does 
not suggest which methodology is more accurate for the 
optically dim, high-redshift sources which dominate in 
the submillimeter, and for which spectroscopic redshifts 
are much rarer. The overall photometric redshift error is 
estimated to be a z ~ 0.06 for both cat alogs. 

CO MBO-17 photometric redshifts (jWolf et all 12004 
120081 are available for the same region of the sky as 
MUSYC. The 5 broad and 12 narrow photometry bands 
from 350 to 930 nm enable accurate SED measurements, 
although the lack of near infrared data limits COMBO- 
17 capabilities to estimate redshifts reliably at z > 1.2. 
The accuracy of these redshifts is a z ~ 0.01 for galax- 
ies with R < 21, 0.02 for galaxies with R ~ 22 and 
0.1 for those with R > 24, which is the sensitivity 
limit o f the survey. Therefo re, following other authors 
(e.g. ILe Floc'h etaH 120(151 ). we only use COMBO-17 
for sourc es with listed R magnitude bri ghter than 24. 

Finally [Rowan- Robinson et alJ (|2008f ) provide a cata- 
log of photometric redshifts for the 24 /xm sources de- 
tected in SWIRE over a larger 4.56 deg 2 region within 
BGS-Wide, but with a shallower depth compared to FI- 
DEL. The optical data (g' r' i') are combined with the 
two shorter IRAC wavelengths, and are used to fit SED 
templates. The accuracy in o~ z is reported to be 0.035 
and is comparable with the above catalogs. 

3. REDSHIFTS OF MIPS SOURCES 

The 9,110 sources in the FIDEL catalog have positional 
accuracies that are significantly smaller than the MIPS 
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PSF (Point Spread Function) at 24 /im. This is because 
the catalog positions are obtained from the more accu- 
rate near infrared (NIR) IRAC catalog. By matching the 
FIDEL and redshift catalogs we identified redshifts for 
~ 72% of the FIDEL sources. We used a search radius of 
1", but using 1.5" or 2" does not substantially change the 
number of identified sources. Given the surface densities 
of sources in the catalogs (MUSYC has one source every 
190arcsec 2 and MUSIC about one every 30arcsec 2 ), we 
expect only a few spurious associations within this search 
radius. 

We identified redshifts from the catalogs in the fol- 
lowing sequence: (i) FIREWORKS (18%); (ii) MUSYC 
(43%); (iii) MUS IC (1%); (iv) COMBO-17 (10 %), and 
(v) the catalog of IRowan-Robinson et all (|2008l ) (0.5%). 
The redshift distribution (solid line in Figure^ peaks at 
z ~ 0.8 and shows a high redshift tail which is particu- 
larly significant in the range 1.5 < z < 2. 

To determine redshifts for the remaining 28% of the 
MIR sources with no direct counterparts in the catalogs 
we used a relation betw een IRAC co l ors an d redshifts 
similar to that derived in iDevlin et al.l (|2009f ) . Redshifts 
from FIREWORKS and MUSYC were matched to IRAC 
sources in the SIMPLE catalog. The resulting -18,000 
objects were binned in a look-up table (LUT) as a func- 
tion of IRAC flux density ratios. The table is shown in 
Figure [3] and exhibits a trend of increasing redshift from 
top- left through top-right to bottom-right. Redshifts for 
galaxies are assigned based on the mean values of sources 
that land in the same bin of the IRAC color-color plane 
with known redshifts. 

A comparison between LUT-retrieved and known red- 
shifts for a number of sources which have not been used 
to compile the table is shown in Figure Although the 
trend has systematic deviations, as well as a large scat- 




FlG. 3. — IRAC color-color plot of SIMPLE sources with red- 
shifts. Pixels are grey-shaded according to redshift with lighter 
shades corresponding to higher redshifts. Overplotted as white 
dots are MIR galaxies in the FIDEL catalog for which no redshift 
information is available. These points cluster in the bottom region 
of the plane, suggesting that they may be a population of higher 
redshift sources. Redshifts were estimated for these objects directly 
from the means of sources with known redshifts that landed in the 
same bins of this color-color plane 




Fig. 4. — Robustness of IRAC redshifts derived from the NIR 
flux density ratios. The top panel compares the IRAC redshifts 
with the catalog of photometric redshifts. Although the relation 
is not linear, its monotonic trend can be used to assign sources to 
redshift bins, provided these are larger than the scatter introduced. 
The redshift binning grid is also shown with the crosses indicating 
the extension of each bin. The distribution of redshifts for sources 
whose redshift lies in the last two bins is shown in the bottom 
panel. The in-bin fraction is 68% and 90%. 

ter, the plot exhibits monotonic behavior, and a reliable 
separation between low and high redshifts. Two plateaus 
are evident, one at 0.6 < z < 1.1 and the other at z > 1.6. 

We initially chose nine redshift bins spanning 0.016 < 
z < 3.5 at logarithmically spaced intervals. The first 
three bins are narrower than the accuracy quoted for the 
redshift catalogs, so we merged them to form a single bin 
4% wide in Az/ (1 + z). We also merged the last two bins 
due to the small numbers of sources that landed within 
them. The six final bins are overplotted in the top panel 
of Figure |4] The redshift distribution for the sources 
selected in the last two bins from the IRAC LUT are 
also shown. Their widths are 68% and 90%, respectively. 
IRAC redshifts are not reliable in the lower redshift bins. 

Almost all of the 28% of the FIDEL sources with no 
redshifts have flux density ratios populating the bottom 
region of the IRAC plane, and their IRAC photomet- 
ric redshift distribution (Figure [2]) indicates that they lie 
preferentially at z > 0.5. This is also confirmed by study- 
ing the distribution of COMBO-17 redshifts for these 
sources. Although unreliable in detail, as they are mostly 
R > 24, the distribution peaks in the same redshift range. 

In summary, we have obtained or estimated redshifts 
for 95% of the sources from the original MIPS catalog. 

4. COMPLETENESS 

The complet eness of the fl u x-lim ited FIDEL catalog is 
investigated in IDevlin et alj (|2009f ). finding it to be 63, 
80 and 96% complete at 20, 40, and 80 /zJy, respectively. 
Ideally, one would like to have a similar completeness 
analysis for each redshift bin used here. Although it is 
tempting to apply the same corrections, in practice one 
should consider that the catalog is based on IRAC posi- 
tions and it would be arbitrary to assume a one-to-one 
relation between sources missing in FIDEL and sources 
missing in the underlying SIMPLE catalog. In this work 
we have therefore ignored any completeness correction. 
However, we do indicate what the expected effect might 
be in the analysis and discussions that follow. 

5. STACKING ANALYSIS 
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TABLE 1 
Stacked Intensities 



A Total Total a (sources Total b (sources with CIB C 

(all sources) with redshift) IRAC-z excluded) 

(im v I v [nWm _2 sr -1 ] 

70 5.4 ±0.3 5.2 ±0.3 4.6 ± 0.3 

250 8.2 ±0.5 8.1 ±0.5 6.1 ± 0.5 10.4 ± 2.3 

350 4.8 ±0.3 4.8 ±0.3 3.3 ± 0.2 5.4 ± 1.6 

500 2.0 ±0.2 2.0 ±0.2 1.2 ±0.1 2.4 ± 0.6 



Note. — Stacked intensities have not been corrected for co mpleteness, and 
the qu oted errors do not include calibration uncertainties. See Marsdcn ct al. 
(2009) for the completeness corrected values. 
a Stacked intensity from every source with a redshift 

b Stacked intensity from every source with a redshift, excluding IRAC- 

estimated redshifts (see text in © 

c FIRAS measured background IFixsen et al.H1998l) 



Stacking analysis using positional information of 
sources selected at different wavelengths is a powerful 
tool for estimating the contribution from a given class of 
objects to the CIB at different wavelengths. In this pa- 
per we have considered the contribution of MIR-selected 
sources to the CIB as a function of redshift at BLAST 
and Spitz er- 70 /xm waveleng ths. Our analysis follows the 
method of lDole et ail ((2006). BLAST maps are whitened 
to suppress the larger scales not relevant to a point- 
source analysis. The filter which makes the noise white 
is estimated from the smoothed two-dimensional power 
spectrum of the map. Scales below ~ 20' are also sup- 
pressed. The whitening filter gain is taken into consider- 
ation when flux densities are measured. The region con- 
sidered conservatively excludes the edges of the FIDEL 
and MUSYC area coverage. This region is described by 
a quadrilateral with corners at the following coordinates 
(J2000): 3 h 31 m 36 s , -28°02'16"; 3 h 31 m 29 s , -27°38'34"; 
3 h 33 m 16 s , -27°35'00"; 3 h 33 m 26 s , -27°58'50". 7,280 FI- 
DEL sources are in this region, compared to the total 
9,110 listed in the full catalog. 

The means of each map in the stacking region are sub- 
tracted, as required by the stacking formulation. We 
then compile a list of postage-stamp maps centered at 
the position of each MIR source. These sub-maps are co- 
added and the total flux density evaluated using aperture 
photometry. Removing the mean locally before stacking 
is formally correct; we verified that the standard tech- 
nique of subtracting the sky estimated in an annulus 
around the aperture gives a compatible result. Whiten- 
ing the maps removes all significant large scale structure 
so that the postage-stamp maps c an be co-ad d ed wi thout 
applying random rotations, as in IDole et all ([2006). Al- 
though the SWIRE 70 /mi map was not whitened, the fil- 
tering applied by the Spitzer reduction pipeline is strong 
enough to remove large scale fluctuations. We verified 
that the retrieved flux densities are the same both by 
stacking postage-stamp maps directly, and by alternately 
rotati ng them by 90° as in the analysis of IDole et all 
((20061) . 

The total retrieved flux density can be overesti- 
mated if the distribu tion of MIR sources is clustered. 
IMarsden et al.l (|2009f ) find this effect to be negligible for 
these catalogs and beam-sizes and it is ignored in the 
subsequent analysis. 

In order to validate the technique and to verify that 



the analysis pipeline does not introduce artifacts, stack- 
ing was performed on noiseless Monte Carlo simulations. 
Sky rea lizations were gener ated at 250 /xm from model 
counts (Laga che et al.1 12003). Source flux densities were 
drawn from the model and redshifts were assigned from 
a uniform distribution spanning < z < 2. Observed 
24 /xm flux densities were assigned to each galaxy using 
Arp 220 as a template to scale from the 250 /xm flux 
densities given the redshifts. Only 24 /xm flux densities 
brighter that 13/iJy were retained in this list. 

We simulated time-stream detector data by scanning 
the sky realizations with the actual BLAST telescope 
pointing solution. We then used the BLAST data re- 
duction pipeline, identical to that for the real data, to 
generate maps. 

Stacking was performed on the simulated sky realiza- 
tions and on the maps obtained from the pipeline, and 
then compared to the total 250 /xm flux density of the 
sources in the 24 /xm mock catalog. The three measure- 
ments were found to agree within 2%, which is well within 
the errors. 

A final test was carried out by stacking the FIDEL 
sources on the 24 /mi Spitzer map itself. Again, we found 
that the retrieved flux density was in good agreement 
with the total flux density in the list. 

The uncertainties are estimated as in [Marsdcn et al.l 
(2009) and a similar check for Gaussianity was per- 
formed. We also verified that stacking against a catalog 
of randomly selected positions results in a value of zero. 

Our stacking region covers an area of only 590arcmin 2 . 
Sampling varianc e can thus play an important role, with 
IDole et all (|2006|) suggesting that it can be as large as a 
factor of 2 (peak-to-peak) . Although we do not presently 
have comparable data in a different field with which to 
test this hypothesis, we make a rough estimate by di- 
viding the present field into 4 separate pieces of roughly 
equal areas. Stacking in each of these regions indepen- 
dently we find an RMS dispersion of ~ 15% in the re- 
trieved flux densities, which is fairly correlated among 
the different BLAST and Spitzer wavebands. Because 
of the 4 sub-fields, the final RMS variation is expected 
to be a factor of 2 smaller than this, although we do 
not explicitly use this estimate in the error bars later in 
our analysis. This is probably only a lower limit on the 
sampling variance as this calculation can not account for 
scales larger than 1 /4 of the map. 
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Table [T] lists the results of this analysis. The 
70 fim retrieved intensit y is compatible with values from 
IDole et all ((200(1 and lDve et all pOOl . At the BLAST 
wavelengt hs these values are within 6% of the ones re- 
p orted bv|Devlin et aLI (|2009f) . and within 4% of those 
in lMarsden et al.1 (|2009l) . 

The stacked intensities for the 5% of sources with no 
rcdshift information is negligible compared to the total 
(< 3%). The table also lists flux densities retrieved us- 
ing only sources with redshifts obtained directly from the 
catalogs (i.e., leaving out the IRAC redshifts). These flux 
densities are 10%, 30%, 30%, and 40% lower at wave- 
lengths 70, 250, 350, and 500 /*m, respectively. This is 
broadly expected if the excluded sources lie at high red- 
shift, such that the longer wavelengths are more seri- 
ously affected. It also provides a consistency check for 
the IRAC rcdshift method discussed earlier. 

Although the stacked flux density as a function of 
wavelength is compatible with the FIRAS detection, 
given its large error and the sampling variance, it is still 
possible that a missing population of faint 24 /im sources 
is needed to resolve the CIB c ompletely. This is men- 
tioned in Ma rsden et al . (2009), and we discuss it fur- 
ther in later sections. In any case, it is worth noting that 
existing estimates of the total resolved CIB are highly 
uncertain, with FIRAS measurements being affected by 
systematics and having a formal error of 25%. 

6. CIB REDSHIFT DISTRIBUTION 

The observed CIB is the total dust reprocessed 
starlight produced in galaxies at all redshifts and lumi- 
nosities. Splitting the FIDEL catalog into redshifts bins 
is an effective way to study the contribution that each 
sub-catalog has to the total intensity. The redshift bins 
used (see §3]) were chosen to be at least a z or larger in 
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Fig. 5. — The cumulative energy distribution is shown as a func- 
tion of 24 fim flux density (left panel) and redshift (right panel) 
for Spitzer at 70 fim, in solid black, and BLAST in blue (250 fim), 
green (350 fJ.ro), and red (500 fim). The 350 fim data are shifted by 
10% to the right in each panel for visual clarity. Each curve shows 
the integrated stacked intensity normalized to the total stacked in- 
tensity at each wavelength as listed in Table [T] Each data point 
is plotted in correspondence of the rightmost edge of each bin. 
Different wavelengths contribute different amounts to the CIB at 
different epochs, or at different 24 (im flux densities, with most of 
the 70 /im background being generated at z < 1 and the 500 fim 
background at higher redshifts. The cumulative plot as a function 
of flux density shows a similar trend at all wavelengths, suggesting 
that brighter 24 fim sources are likely to be at lower redshifts. The 
cumulative energy distribution for Spitzer at 24 fim is also shown 
for reference (dashed line in the right panel). 



Az/(1 + z), and the width of the last two bins ensures 
that the errors in the IRAC redshifts do not result in 
sources being placed in the wrong bin (Figure [J). The 
mean stacked intensity in each bin is listed in Table [2] 

The cumulative distribution of stacked intensities is 
shown in Figure[5]as a function of both 24 /im flux density 
and redshift, and it is normalized to the total intensity 
listed in Table [TJ We prefer to show this normalization 
rather than normalizing to the FIRAS measurements be- 
cause: i) there is no 70 /im FIRAS detection and the 
most accurate detections of the CIB at this wavelength 
are lower limits; ii) the FIRAS detections in the BLAST 
bands have large uncertainties; and iii) the BLAST sur- 
vey is affected by sampling variance, making it difficult 
to compare with an all sky survey. 

A similar qualitative behavior is seen in both panels 
of Figure El with the stacking measurements at shorter 
wavelengths dominating the brighter 24 /im flux densities 
and the lower redshifts. By redshift 1.1, about 75% of 
the CIB is generated at 70 /im, 55% at 250 /im, 45% at 
35 ^m, and 40% at 50 ^m. Similar results are found 
by lMarsden et all (2009) who estimate these fractions to 
be 50, 45, and 40% at 250, 350, and 500 /im, respec- 
tively (although z = 1.2 is the dividing redshift in their 
study) . The relatively small discrepancy at 250 /im can 
be accounted for by the different way redshifts are as- 
signed. In this work, most redshifts have b een retrieved 
from catalogs, while iMarsden et al.l (|2009[ ) only uses an 
approximate statistical division into two redshift bins us- 
i ng IRAC colors. 

iWang et al] (|2006l) stacked F-band and IRAC selected 
sources in the 850 /im map of the GOODS-N field. They 
find that about 70% of the total retrieved flux originates 
at redshifts z < 1.5. This is about 30% of the FIRAS 
measurement, suggesting that a large fraction of the CIB 
at this wavelength m ay be generated at higher redshifts. 
iSerieant et al.l (|2008l ) reach similar conclusions stacking 
MIR and NIR selected sources. 

It is worth mentioning that the fractions of the resolved 
background depend on the normalization adopted. If 
FIRAS is used instead, all of these fractions change by as 
much as 20% at the BLAST wavelengths (not accounting 
for the quoted measurement uncertainties). 

Compared to 70 /im, the 500 /im cumulative distribu- 
tion shows no evidence of convergence in the highest red- 
shift bin, suggesting that a population of faint 24 /mi 
sources is missing in order to fully resolve the CIB at 
these longer wavelengths. 

7. AVERAGE PHYSICAL PARAMETERS 

The average flux densities at 70, 250, 350, and 500 fim 
from stacking (total stack divided by the number of con- 
tributing sources) were used to fit SEDs in each redshift 
interval to estimate the average FIR luminosities (from 
8 to 1000 /im in the rest frame) of the MIR sources as a 
function of redshift. 

Since the emission mechanism at these wavelengths is 
thermal, the natural choice of SED is a modified grey- 
body with some emissivity law: A B{y, T), where 
B(u,T) is the blackbody spectrum and A its amplitude. 
Such a model assumes a single dust temperature, but 
the reality is far more complex. Temperature distribu- 
tions are observed both within single galaxies and among 
different sources. As a consequence, the resulting spec- 
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Fig. 6. — SED fitting of the average flux densities measured in different redshift bins. The points with error bars are from BLAST 
(color-corrected 250, 350, and 500 ^m), and Spitzer 70/xm. The solid lines are the best fits at each redshift, with the 68% confidence levels 
indicated by the dashed lines. The fit accounts for the finite BLAST bandwidths and for the correlated calibration uncertainties. The 
model template is a modified greybody with an emissivity law j3 = 1.5 and a power law v~ a replacing the Wien part of the spectrum. Two 
templates were used: one with a = 2.8, describing Arp 220's spectrum, and one with a = 1.8 which is better matched to M82. Both models 
are compatible with the data, although for clarity only the a = 2.8 fit is shown here. The 24 (im average flux densities of the sources in 
the catalog are also shown in each panel for completeness, but are not used in the fits because these fluxes come from regions which are 
not in thermal equilibrium. 



trum is better described by a power law than a Wien 
exponential decay at wavelengths shorter than the peak 
in the SED. Indeed, for our data a modified greybody 
alone gives a poorer fit to the combined Spitzer 70 /im 
and BLAST data points as redshift increases; this is 
because the 70 ^m channel samples shorter rest-frame 
wavelengths. 

Instead of a blackbody curve modified purely on the 
Rayleigh Jeans side, we instead used a modified grey- 
body with a fixed emissivity index (3 — 1.5, together with 
a power law decay v~ a at short wavelengths to prevent 
the high-frequency SED from falling exponentially. The 
exponent a is chosen by fitting the SEDs of two often 
used galaxies: the Ultra Luminous IR Galaxy (ULIRG, 
L > 1O 12 L , a = 2.8) Arp 220, and the Luminous IR 
galaxy (LIRG, L > 10 n L Q , a = 1.8) M82. The former 
is a merging system, and the latter is a starburst proba- 
bly triggered by tidal encounters with the nearby spiral 
M81. Our fitting procedures estimate the amplitude of 
the template and the temperature, keeping a and [3 fixed. 
Both SEDs fit the four data points at each redshift, as 
shown in Figure[6]for a = 2.8; the best-fit quantities are 
listed in Table H 

The estimate of FIR luminosity (the integral from 8 to 
1000 /jm of the SED in the rest-frame) is a weak function 
of the template (Figure [7|) , emphasizing that the model 
need only provide a reasonable interpolation of the data 



around the peak of the rest-frame FIR emission. The re- 
trieved luminosities increase steeply with redshift. The 
top panel in Figure [7] plots the mean FIR luminosity as 
a function of the redsh ift. It is in good agreem ent with a 
similar result found bv lLe Floc'h et all (|2005l ). who stud- 
ied the evolution of FIR luminosities for the same pop- 
ulation of 24 /jm sources, although with a higher MIR 
flux density limit of 80 /iJy. The figure provides an esti- 
mate of the selection effect arising from the flux density 
limit of the FIDEL catalog, which we assume is 20 /iJy, 
the 63% completeness limit. The calculation uses tem- 
plate models for "no rmal" and "starburst" galaxies from 
iLagache et al.1 (|2003l ). co-added after first normalizing to 
the same FIR l uminosity. 

As shown in iKennicuttl £1998), the FIR emission is 
tightly correlated with star formation activity (Star For- 
mation Rate, or SFR), under the assumption that star 
formation is completely optically obscured by dust and 
that the fraction of infrared l uminosity d omina ted by 
AGN heating is neglig ible 13 . iKennicuttJ l)1998[ ) shows 
that the relation between FIR luminosity and dust ob- 
scured SFR is then: 



SFRpvloyr- 1 ] = 1.728 x 10~ 10 L [L ]. 

13 IDevlin et all H2009I ) show that the AGN contribution to the 
CIB at BLAST wavelengths is about 7%, and it is small enough to 
be ignored. 
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Fig. 7. — The mean rest-frame luminosities, and temperatures of 
the MIR FIDEL sources are plotted against redshift. In the top 
panel, the mean FIR luminosities are estimated from the a = 2.8 
(black points) and a = 1.8 (grey, dashed points) SED fits. There 
is little difference between the two, confirming that the estimated 
luminosities are a weak function of the template used. The grey 
dashed line provides an estimate of selection arising from a 24 (im 
flux density limit of 20 fijy (63% completeness limit of the cata- 
log). The bottom panel shows rest-frame temperatures for both 
templates. A difference of ~ 5 K is observed between the two, with 
temperatures apparently increasing with redshift. 

We use this relationship to calibrate the SFR axes of the 
luminosity plots. 

From Figure [7J we observe that at low redshifts (z < 
0.6) the MIR population is dominated by moderate IR 
emitters (L < 6 x 10 10 L Q ). As the redshift increases, the 
volume probed also increases and LIRGs become more 
common, which is consistent with many spheroids being 
formed at about this time [z ~ 1). This is expected, 
as SCUBA surveys at 850 /im have revealed the exis- 
tence of a large number of ULIRG-type objects, residing 
at typically higher re dshifts with a median of z ~ 2.5 
(jChapman et al.ll2003fl . 

Stacking analyses cannot be used to study individ- 
ual sources, only their average properties. Therefore 
we emphasize once more that the quantities estimated 
are the average values of the population of MIR sources 
resolving most, if not all, of the CIB detected by FI- 
RAS. For instance, the last bin in Figure [JJ probes a red- 
shift range similar to the SCUBA galaxies. However, 
the mean luminosity is slightly lower than that typically 
quoted for SCUBA galaxies, suggesting either that some 
of the most luminous objects (and potentially most opti- 
cally obscured) are missing from the bin, or that SCUBA 
galaxies have luminosities that are above average at those 
redshifts (which is surely true for the extreme objects 
that have been studied so far at 850 /im) . 

Our template fitting returns rest-frame temperatures, 
but these are more model-dependent than the luminosi- 
ties. As shown in Figure [JJ the difference in the tem- 
peratures determined for the two templates is about 5 K 
across the whole redshift range. 

We observe an increase in mean temperature with in- 
creasing redshift. This trend can be explained by the 



10 1 



10 9 r 



10 r 
10 7 r 



10 6 r ? 



10 c 



10 b 



t 



10 1 



It 



10 1 



10 1 



Luminosity [L ] 



Fig. 8. — The mean dust masses of the MIR FIDEL sources 
are plotted against the mean FIR luminosities of Figure [7] These 
masses are calculated using the temperatures estimated from the 
a = 2.8 (black points), and a = 1.8 (grey, dashed points) SED fits. 
Redshifts are not indicated explicitly, but increase monotonically 
from left to right. 



Fig. 9. — Comoving dust mass density vs. redshift. There is only 
moderate evidence of evolution with redshift. 

correlation observe d between temperat ure and luminos- 
ity as discussed by IDunne et al.1 (|2000r i. Detailed mod- 
eling beyond the scope of this paper would be required 
to determine whether selection or cosmic evolution also 
play a role. 

The average mass in dust can be estim ated from the 
observed flux densities (Hildebrand 1983), S v (vq): 



Ma = 



dm 

1 + Z K (v em ) B v (Vem, T) 



where Z?l is the luminosity distance. The mass- 
absorption coefficient, k, is evaluated at 1 mm (Ki mm = 
0.1m 2 kg" 1 , IHughed [T996). and has a v" dependence. 
Any dust mass estimate has to be treated with cau- 
tion becaus e: i) Ki mm is un certain within an order of 
magnitude (|Blain et al.ll2002T ); and ii) the temperature is 
largely uncertain, and model dependent. We use the flux 
density at a rest-frame wavelength of 500 /im, and tem- 
peratures obtained from the two model SED fits. The 
average dust masses are plotted against the average lu- 
minosity in Figure [51 Given the monotonic relation be- 
tween FIR luminosities and redshifts, the leftmost and 
rightmost points in the figure are the lowest and highest 
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Redshift 

Fig. 10. — The evolution of the total FIR luminosity density up 
to z = 3.5 from MIR-selected sources is shown as blac k filled cir- 
cles. T he other points are optical-UV detections fromlLillv et al.l 
(1996) (up-pointing solid triangles) and Stcidcl ct al. (1999) (left- 
pointing solid triangles). These measurements need to be corrected 
upward for the effects of dust extinction before they can be com- 
pared to detections in the FIR, and the corresponding symbols with 
dashed error bars are th e extinction-correcte d values. The open 
circle is the lower limit of Hughes ct al. ( 1998) derived for SCUBA 
galaxies. SFR is shown on the left axis and shows good agree- 
ment with the optical extinction corrected estimates, confirming 
that most of the star-forming activity in the Universe is obscured 
by dust. The light grey filled circles, with dashed error bars, are 
calculated from the same MIR population, from which the sources 
with redshifts estimated from IRAC colors have been removed, and 
are shifted by 10% toward increasing z for visual clarity. A luminos- 
ity function model (not a fit) is plotted as a solid line, and it is also 
evaluated at a 24 (im flux density limit of 20 fijy to indicate the ef- 
fect of the FIDEL catalog flux density limit (da shed line) . The data 
also s how good agreement with the model of [Hopkins & Beacoml 
(2006) (dash— dot line). The last redshift bin falls low compared 
with the SCUBA lower limit (and the extinction-corrected point 
at similar redshift), and suggests that a fraction of SCUBA galax- 
ies might be a population of faint 24 (im sources missing from the 
FIDEL sample. 

redshift bins, respectively. The dust mass shows a steady 
increase over two orders of magnitude across the probed 
range of redshifts, and FIR luminosities. 

We observe that the estimates of temperature and dust 
mass from the fit with the a = 1.8 template (M82) 
have larger errors compared to those obtained from the 
a = 2.8 template (Arp 220). This is due to the flat- 
ter spectrum close to the emission peak of the a = 1.8 
template. 

8. LUMINOSITY AND STAR FORMATION HISTORIES 

Mean rest-frame FIR luminosities and dust masses can 
be converted into comoving-volume densities. The con- 
version factor applied to each redshift bin is the ratio of 
the number of sources in the bin to the comoving-volume 
probed by the bin, which, for a Euclidean geometry, we 
calculate as 



tack 



gj(fhO d 3 M 

(1 + Z hi ) 3 (1 + 2lo) 3 



Here, ^ s tack is the solid angle of the stacking region, z\ 
and Zhi are the edges of each redshift interval, and z is 
the bin center. 

Interestingly, the comoving dust mass densities plotted 
in Figure[H]show little variation over the redshifts probed. 
This suggests that the amount of comoving dust associ- 
ated with the FIDEL sources is fairly constant as galaxies 



are built up. iDunne et al.l (2003) found that there was 
much more dust at z > 1, although their survey is only 
sensitive to galaxies with very large dust masses, whereas 
our stacking analysis is sensitive to galaxies with a wider 
range of dust masses. 

The comoving FIR luminosity densities are shown in 
Figure [TU] and are listed in Table O The left vertical 
axis is converted into star formation rate density. We 
used the template with a = 2.8 to calculate the lumi- 
nosities. Also shown is the effect of excluding sources 
for which redshifts were obtained from the IRAC color- 
color plane. This potential incompleteness mostly affects 
the higher redshift bin. We note once again that the er- 
ror bars include only the statistical uncertainties arising 
from the stacking and calibration and do not include the 
effect of sampling variance or incompleteness. This plot 
is a complete representation of the energetics in the CIB 
generated by MIR sources at the FIDEL flux density 
limit, as 95% of these galaxies have been included and 
the remaining 5% do not contribute substantially to the 
background. The luminosity density 

plot shows the luminosity integrated over the comoving 
volume. A naive model is overplotted ( not fitted) and is 
derived from the I Saunders et al.l (|1990f ) luminosity func- 
tion (LF) form for IRAS galaxies, 



exp 



1 



2^r 



log 2 1 + — 



T * 



The coefficient s ado pted for this model are from 
iLe Floc'h et~aTl (|2005l ) for a fit to the local FIR LF. The 
quantity plotted is then 



dL 



L $IR 



L 

W) 



dlogL. 



L m ; n is set to zero in one case (dashed curve) to show 
the total FIR luminosity density predicted by this model. 
The solid curve was calculated for L m - ln estimated assum- 
ing a 24 /jm flux density of 20 /xJy (63% completeness) 
and our fitting template, and shows the effect of the flux 
density limit in the FIDEL catalog (higher flux density 
limits correspond to larger values of L min ). The evolu- 
tion was chosen to be f(z) = (1 + z) 3 2 up to z = 1.2, 
and constant at higher redshifts. This naive modeling 
shows good agreement with the data and qualitatively 
demonstrates the effect of the flux density limit of the 
catalog. As expected, the only substantial change is to 
the highest redshift point. 

The FIR luminosity densities estimated here show 
good agreeme nt with similar re sults obtained from opti- 
cal detections (jLillv et al.|[l99 6) at z < 1, after a correc- 
tion for dust extinction has been applied by the authors. 
Optical measurements are sensitive to the fraction of 
photons not absorbed by dust. Therefore dust extinction 
needs to be taken into consideration in order to estimate 
the total star fromation activity, and to compare with 
measurements obtained in th e FIR. In Figure [TO] we also 
plot the empirical model of Hopkins & Bcacom (2006) 
obtained from a parametric fit to SFR densities measured 
at wavelengths spanning from the radio to the UV, which 
shows good agreement with our measurements. 

As the MIR sources resolve most, if not all, of the CIB, 
at its emission peak, this is the first direct measurement 
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TABLE 2 

Stacking Analysis Results 







vl v (70 fim) vl v (250 /im) ul v (350 fim) vl v (500 fim) 




Temperature a 


Dust mass a 


LlR density Sources in bin 






nWm _2 sr _1 


10 9 L o 


K 


1O 7 M 


10 8 L o Mpc" 3 



0.016 0.098 0.53 ±0.04 0.31 ± 0.07 0.11 ±0.04 0.05 ± 0.02 1.2±%\ 23.2±o_| [21.7^^;!] 1.8±f% [2.3^™] 1.19±g;J| 115 E 



0.098 0.177 0.49 ±0.04 0.36 ± 0.08 0.13 ±0.04 0.06 ± 0.03 5.3t ;| 24.2±?;§ pO.ajl^a] 8 - 4 -5°i [H^ljfg 2 ] l- 6 9±o.i3 169 

0.177 0.322 0.54 ±0.06 0.20 ±0.11 0.13 ±0.06 0.06 ± 0.03 9.71JJ 29.5±g 3 [26.2t|;f] 7.0jl 2 , ;| I 9 - 4 !^] i-Hto.is 323 

0.322 0.585 1.31 ±0.11 1.19 ±0.19 0.45 ±0.10 0.19 ±0.06 41.3±!;} 29.8±o_| [26.7±5".i] 32 - 9 ±2i'.9 [ 49 - 2 -3o] 2 - 96 to.i3 949 

0.585 1.062 1.26 ±0.16 2.53 ± 0.30 1.38 ±0.16 0.42 ± 0.09 103±g° 30-Q±S;I [24.3±2;J] 128tgo 4 [236±?gg] 4.54± ;|2 2247 

1.062 3.500 1.06 ±0.19 3.45 ± 0.35 2.51 ±0.18 1.22 ±0.11 899±gg 36A±f° ^.l 4 ;^] 936l^3 8 [1363175™] 6 - 39 ±o.53 3115 



Note. — FIR luminosities, Temperatures and dust masses are the mean values per source in each redshift interval. The FIR luminosity densities are calculated 
from the mean FIR luminosities as described in Section [8] 
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of its energetics and of the dust obscured star forma- 
tion history of the Universe. Although the agreement 
with previous measurements is remarkably good, previ- 
ous studies have had to rely on modeling of dust obscu- 
ration and SEDs to evaluate FIR luminosities. With the 
combination of BLAST and Spitzer, we have been able 
to directly probe this quantity with very few assumptions 
about the underlying physical details. 

Comparing our measurements with the optical mea- 
surements (the triangle symbols with solid error bars 
in Figure [TO)) , we find that the fraction of light from 
young stars hidden by dust is about 70% in the range 
0.1<z<l. 

The highest redshift bin, which coincides with the 
redshifts at which SCUBA galaxies are typically found, 
is the most affected by the FIDEL flux density limit. 
Comparing with the lower limit o n the S FR density for 
SCUBA galaxies of Hughe s et al.l ([1998), it appears as 
if some fraction of this population is missing from the 
FIDEL catalog, which does not contain large fractions of 
ULIRG-type objects, as discussed in the previous section. 
It is possible that a significant fraction of SCUBA galax- 
ies have fainter 24 /im flux densities than the threshold we 
used for the FIDEL catalog, although we note that many 
SCUBA sources ha ve had 24 /zm counterparts identifie d 
in past surveys (e.g. lPope et al.l l2006; Ivi son et al.ll2007l) . 
However, here we are talking about the sources domi- 
nating the background at 850 /im, which are fainter than 
the typical sources detected by SCUBA. This is indicated 
by our simple luminosity function model, which suggests 
that a substantial number of low-luminosity sources are 
missing from the catalog at high redshifts. It also may 
be possible that some of the high-redshift sources in our 
catalog have simply been mis-identified as lower-redshift 
objects, but due to their relatively small numbers, have 
little effect on the lower-redshift bins. 

Either way, there is evidence that this study misses 
a portion of the SFR history at the highest redshifts 
(z > 1), during which it is believed the most massive 
galaxies were forming throug h powerful mergers. T his is 
consistent with the result of Marsdcn et al. (2009) who 
find the contribution of FIDEL sources to the CIB at 
850 /im to be significantly low compared to the FIRAS 
measurements. They also find that their estimate of the 
CIB generated at 850 /im b y sources with z > 1. 2 is low 
compared to the model of Valiantc et al.l (|2009f) which 
otherwise fits the CIB at the BLAST wavelengths. 

9. CONCLUSIONS 

We have studied the contribution of MIR galaxies se- 
lected at 24 /im to the CIB and found that they resolve 
most, if not all, of the radiation detected by the all- 
sky FIRAS and DIRBE surveys at the peak of the far- 
infrared background. We assigned redshifts to more than 
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70% of this population from existing redshift catalogs, 
and developed a technique which uses the NIR IRAC 
colors to statistically assign redshifts to the remainder of 
the catalog. Using a stacking analysis to overcome the 
confusion noise arising from the finite instrumental reso- 
lution, we study the composition of the CIB as a function 
of redshift. We find that 60% of the CIB at 500 /im is 
generated at redshifts z > 1.1, while the 70 /im back- 
ground has a more recent origin, with 75% of it being 
generated at redshifts z < 1.1. 

By fitting SEDs to the mean flux densities at each red- 
shift, we have shown that, on average, this population 
is consistent with being predominantly lower-luminosity 
LIRGs rather than the ULIRGs detected in SCUBA sur- 
veys. It appears that a significant fraction of the sub- 
millimeter galaxies known to exist at redshifts z > 2 
are missing from our analysis, and we conjecture that 
they are simply below the 24 /im flux density limit in 
our source catalog, or perhaps in some cases have had 
their redshifts mis-identified. However, while such ob- 
jects are important at high redshifts, their contribution 
to the peak of the CIB is small compared with the bulk 
of the MIR-selected galaxies at lower redshifts. 

The evolution in the total comoving FIR luminosity 
density can be evaluated using MIPS and BLAST data 
directly, as these wavelengths span the rest-frame FIR 
peak emission. This evolution can then be directly con- 
verted into the star formation rate history. We find good 
agreement with existing studies conducted at optical-UV 
wavelengths which have had to rely on many more as- 
sumptions about physical source properties. 

Our results confirm that dust obscuration dominates 
the history of star formation of the Universe, with star 
formation rates that are about three times larger than in 
the optical-UV, in the redshift range 0.1 < z < 1. 

These results can be used to constrain models of lu- 
minosity densities and galaxy formation as they directly 
probe the FIR energetics of the Universe up to redshift 
z ~ 3. 
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